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ABSTRACT 



CJ ■ We present a fully covariant and gauge-invariant calculation of the evolution 

O , of anisotropics in the Cosmic Microwave Background (CMB) radiation. We use 

the physically appealing covariant approach to cosmological perturbations, which 
^ . ensures that all variables are gauge-invariant and have a clear physical interpre- 

tation. We derive the complete set of frame-independent linearised equations 



O ■ describing the (Boltzmann) evolution of anisotropy and inhomogeneity in an al- 

^ . most Friedmann- Robert son- Walker (FRW) Cold Dark Matter (CDM) universe. 

O , These equations include the contributions of scalar, vector and tensor modes in 

Q\ . a unified manner. Frame-independent equations for scalar and tensor perturba- 

tions, which are valid for any value of the background curvature, are obtained 
straightforwardly from the complete set of equations. We discuss the scalar equa- 
tions in detail, including the integral solution and relation with the line of sight 



c/3 , approach, analytic solutions in the early radiation dominated era, and the numer- 

^ ' ical solution in the standard CDM model. Our results confirm those obtained by 



other groups, who have worked carefully with non-covariant methods in specific 
gauges, but are derived here in a completely transparent fashion. 

Subject headings: cosmic microwave background — cosmology: theory — gravi- 
tation — large-scale structure of universe 



1. Introduction 



The cosmic microwave background radiation (CMB) occupies a central role in modern 
cosmology. It provides us with a unique record of conditions along our past lightcone back 
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to the epoch of decouphng (last scattering), when the optical depth to Thomson scattering 
rises suddenly due to Hydrogen recombination. Accurate observations of the CMB anisotropy 
should allow us to distinguish between models of structure formation and, in the case of non- 
seeded models, to infer the spectrum of initial perturbations in the early universe. Essential 
to this programme is the accurate and reliable calculation of the anisotropy predicted in 
viable cosmological models. 

Such calculations have a long history, beginning with Sachs & Wolfe (1967) who investi- 
gated the anisotropy on large scales (;^ 1°) by calculating the redshift back to last scattering 
along null geodesies in a perturbed universe. On smaller angular scales one must address the 
detailed local processes occurring in the electron/baryon plasma prior to recombination, and 
the effects of non-instantaneous last scattering. These processes, which give rise to a wealth 
of structure in the CMB power spectrum on intermediate scales and damping on small scales 
(see, for example. Silk (1967, 1968)), are best addressed by following the photon distribution 
function directly from an early epoch in the history of the universe to the current point 
of observation. This requires a numerical integration of the Boltzmann equation, and has 
been carried out by many groups, of which Peebles & Yu (1970), Bond & Efstathiou (1984, 
1987), Hu & Sugiyama (1995), Ma & Bertshcinger (1995), Seljak & Zaldarriaga (1996) is a 
representative sample. 

The calculation of CMB anisotropics is simple in principle, but in reality is plagued 
with subtle gauge issues (Stoeger, Ellis, & Schmidt 1991; Stoeger et al. 1995; Challinor & 
Lasenby 1998). These problems arise because of the gauge- freedom in specifying a map $ 
between the real universe (denoted by S) and the unperturbed background model(denoted by 
S) [(Ellis fc Bruni 1989 )| , which is usually taken to be a Friedmann-Robertson- Walker (FRW) 
universe. The map $ identifies points in the real universe with points in the background 
model, thus defining the perturbation in any quantity of interest. For example, for the 
density p as measured by some physically defined observer, the perturbation at x G S* is 
defined to be 6p{x) = p{x) — p(x), where p is the equivalent density in the background 
model, and x maps to x under $. The map $ is usually (partially) specified by imposing 
coordinate conditions in S and S. Any residual freedom in the map $ after the imposition 
of the coordinate conditions (gauge-fixing) gives rise to the following gauge problems: (i) the 
map cannot be reconstructed from observations in S alone, so that quantities such as the 
density perturbation, which depend on the specific map $, are necessarily not observable; 
(ii) if the residual gauge freedom allows points in S to be mapped to physically inequivalent 
points in S in the limit that S = S, then unphysical gauge mode solutions to the linearised 
perturbation equations will exist. 

There are several ways to deal with the gauge problems described above. In the earliest 
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approach (Lifshitz 1946) , one retains the residual gauge freedom (in the synchronous gauge) 



but keeps track of it so that gauge mode solutions can be eliminated. Furthermore, the final 
results of such a calculation must be expressed in terms of the physically relevant, gauge- 
invariant quantities. Although there is nothing fundamentally wrong with this approach if 
carried out correctly, it suffers from a long history littered with confusion and errors. The 
need to express results in terms of gauge-invariant variables suggests that it might be benefi- 
cial to employ such variables all along as the dynamical degrees of freedom in the calculation. 
A further advantage of such an approach is that gauge modes are automatically eliminated 
from the perturbation equations when expressed in terms of gauge-invariant variables. This 
is the approach adopted by Bardeen (1980), who showed how to construct gauge-invariant 
variables for scalar, vector and tensor modes in linearised perturbation theory, by taking 
suitable linear combinations of the gauge- dependent perturbations in the metric and matter 
variables. This approach has been used in several calculations of the CMB anisotropy (see, 
for example, Abbott & Schaefer (1986) and Panek (1986)). However, the Bardeen variables 
are not entirely satisfactory. The approach is inherently linear, so that the variables are 
only defined for small departures from FRW symmetry. Furthermore, the approach assumes 
a non-local decomposition of the perturbations into scalar, vector and tensor modes at the 
outset, each of which is then treated independently. As a result, the Bardeen variables are 
only gauge-invariant for the restricted class of gauge-transformations that respect the scalar, 
vector and tensor splitting. Finally, although the Bardeen variables are gauge-invariant, 
they are not physically transparent, in that, in a general gauge, they do not characterise the 
perturbations in a manner that is amenable to simple physical interpretation. 

An alternative scheme for the gauge-invariant treatment of cosmological perturbations 
was given by Ellis & Bruni (1989) (see also, Ellis, Hwang, & Bruni (1989)) who built upon 
earlier work by Hawking (1966). In this approach, which is derived from the covariant 
approach to cosmology/hydrodynamics of Ehlers and Ellis (Ehlers 1993; Ellis 1971), the 
perturbations are described by gauge-invariant variables that are covariantly defined in the 
real universe. This ensures that the variables have simple physical interpretations in terms 
of the inhomogeneity and anisotropy of the universe. Since the definition of the covariant 
variables does not assume any linearisation, exact equations can be found for their evolution, 
which can then be linearised around the chosen background model. Furthermore, the co- 
variant approach does not employ the non-local decomposition into scalar, vector or tensor 
modes, at a fundamental level. If required, the decomposition can be performed at a late 
stage in the calculation to aid solving the equations. Even if one denies that working with 
gauge-invariant variables is a significant advantage, the key advantage of the covariant ap- 
proach, however, is that one is able to work exclusively with physically relevant quantities, 
satisfying equations that make manifest their physical consequences. 
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The covariant and gauge-invariant approach has aheady been apphed to the hne of 
sight calculation of CMB anisotropics under the instantaneous recombination approximation 
(Dunsby 1997; Challinor & Lasenby 1998), and has been used to obtain model- independent 
limits on the inhomogeneity and anisotropy from measurements of the CMB anisotropy on 
large scales [(Maartens, Ellis, fc Stoeger 1995 )| . In this paper, we extend the methodology de- 



veloped in these earlier papers, to give a full kinetic theory calculation of CMB anisotropics 
valid on all angular scales. Our motivation for reconsidering this problem is two-fold. Firstly, 
it is our belief that the covariant and gauge-invariant description of cosmological perturba- 
tions provides a powerful set of tools for the formulation of the basic perturbation equations, 
and their subsequent interpretation, which are superior to the techniques usually employed 
in such calculations for the reasons discussed above. In particular, by applying covariant 
methods for the problem of the generation of CMB anisotropics, we can expect the same 
advantages of physical clarity and unification that have already been demonstrated in other 
areas, (Ellis et al. 1989; Bruni, Ellis, & Dunsby 1992; Dunsby, Bruni, & Ellis 1992; Dunsby, 
Bassett, & Ellis 1996; Tsagas & Barrow 1997). The approach described here brings the 
underlying physics to the fore, and can only help to consolidate our rapidly growing un- 
derstanding of the physics of CMB anisotropics. Furthermore, although we only consider 
the linearised calculation here, the extension of these methods to the full non-linear case is 
quite straightforward (Maartens, Gebbie, & Ellis 1998). Our second motivation is to perform 
an independent verification of the results of other groups (for example. Ma & Bertschinger 
(1995)), with a methodology that is free from any of the gauge ambiguities that have caused 
problems and confusion in the past. Given the potential impact on cosmology of the next 
generation of CMB data, we believe that the above comments provide ample justification 
for reconsidering this problem. 

For definiteness we consider the cold dark matter (CDM) model, although the methods 
we describe are straightforward to extend to other models. We have endeavoured to make 
this paper reasonably self-contained, so we begin with a brief overview of the covariant ap- 
proach to cosmology and define the key variables we use to characterise the perturbations in 
Section ^ We then go on to present a complete set of frame-independent equations describ- 
ing the evolution of the matter components and radiation in Section § in an almost-FRW 
universe (with arbitrary spatial curvature). These equations, which employ only covariantly 
defined, gauge-invariant variables, are independent of any harmonic analysis; they describe 
scalar, vector and tensor perturbations in a unified manner. Many of the equations have sim- 
ple Newtonian analogues, and their physical consequences are far more transparent than the 
equations that underlie the metric-based approaches. Equations pertinent to scalar modes, 
see Section ^ and tensor modes, see Section |^, can be obtained from the full set of equations 
with very little effort, and are useful at this late stage in the calculation as an aid to solving 
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the linearised equations. A significant feature of this approach is that a covariant angular 
decomposition of the distribution functions is made early on in the calculation, before any 
splitting into scalar, vector and tensor modes. This allows scalar, vector and tensor modes 
to be treated in a more unified manner. In particular, the azimuthal dependence of the 
moments of the distribution functions does not have to be put in by hand (after inspection 
of the azimuthal dependence of the other terms in the Boltzmann equation), as happens 
in most metric-based calculations. This is particularly significant for tensor modes where 
the required azimuthal dependence is non-trivial and is different for the two polarisations of 
gravitational waves. We consider the equations for scalar modes in considerable detail. We 
present the integral solution of the Boltzmann multipole equations in a = almost-FRW 
universe, and discuss the relation between line of sight methods (which employ lightlike 
integrations along the lightcone) and the Boltzmann multipole approach (where a timelike 
integration is performed). We derive analytic solutions for scalar modes in the early radia- 
tion dominated universe, that are used as initial conditions for the numerical solution of the 
scalar equations, the results of which we describe in Section p. In Section |^ we give a brief 
discussion of the tensor equations in the covariant approach. The covariant angular decom- 
position naturally gives rise to a set of variables that describe the temperature anisotropy in 
a more direct manner than in the conventional metric-based approaches. This is particularly 
apparent for tensor perturbations, where the CMB power spectrum at a given multipole / 
is determined by the / — 2, / and (/ + 2)-th moments of the conventional decomposition of 
the photon distribution function, which obscures the physical interpretation of these mo- 
ments. Finally, we end with our conclusions in Section ^ Ultimately, our results confirm 
those of other groups (for example. Ma & Bertschinger (1995)) who have performed similar 
calculations by working carefully in specific gauges, but are obtained here with a unified 
methodology that is more physically transparent and less prone to lead to confusion over 
subtle gauge effects. 

We employ standard general relativity and use a (H ) metric signature. Our 

conventions for the Riemann and Ricci tensors are fixed by [Va, Vb\u'^ = —TZabd'^u'^, and 
TZab = Tlacb'^- Round brackets around indices denote symmetrisation on the indices enclosed, 
and square brackets denote antisymmetrisation. We use units with c = G = 1 throughout, 
and a unit of distance of Mpc for numerical work. 

2. The Covariant Approach to Cosmology 

In this section, we summarise the covariant approach to cosmology (Ehlers 1993; Ellis 
1971; Hawking 1966) and the gauge- invariant perturbation theory of Ellis & Bruni (1989) 
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which is derived from it. We begin by choosing a velocity field u"- which satisfies the following 
criterion: the velocity must he physically defined such that it reduces to the four velocity of 
the fundamental observers in the Friedmann-Robertson-Walker (FRW) limit. This restriction 
on -u" is essential to ensure gauge-invariance of the Ellis & Bruni perturbation theory. Note 
that, in a general perturbed spacetime, there is no unique choice for -u". Acceptable choices 
for u"" include the four velocity of a given matter component, and the timelike eigenvector of 
the stress-energy tensor. In the covariant approach to perturbations in cosmology, covariant 
variables are introduced that describe the inhomogeneity and anisotropy of the universe. 
These variables employ the velocity field u"- in their definition, and so, in a given spacetime, 
their values depend on how we choose (the exact transformation laws are given in Maartens 
et al. (1998)). For a given choice of u"-, the covariant variables, defined below, describe the 
results of observations made by observers comoving with the velocity m", and their frame- 
dependence reflects the fact that the observations depend on the velocity of the observer. 
It might be thought that the freedom in the choice of velocity would introduce similar 
ambiguities as the choice of map $ does in conventional approaches. However, this is not 
the case because of the restriction on u"- that we emphasised above. It is certainly true 
that, with a suitable choice of u"-, we can eliminate some aspect of the inhomogeneity and 
isotropy observed. For example, we can always choose u"' so that the CMB dipole vanishes. 
However, in a given spacetime, the covariant variables cannot be forced to take arbitrary 
values through some particular choice of u"". In particular, if, for some timelike velocity 
field (not necessarily restricted to satisfy the criterion emphasised above), all of the gauge- 
invariant variables defined below vanish identically, then the universe is necessarily FRW. 
This gives a covariant condition which characterises the FRW limit, but note that, if u"' is 
unrestricted, we could have the situation where the universe is FRW, but we are not viewing 
it from the perspective of the fundamental observers and so some of the variables would 
not vanish. (This is similar to the presence of gauge mode solutions in the metric-based 
approach.) However, if we ensure that u"" is defined physically, so that in the FRW limit 
it necessarily reduces to the velocity of the fundamental observers, this situation cannot 
arise, and the variables used to characterise the anisotropy and inhomogeneity are genuinely 
gauge-invariant . 

We refer to the choice of velocity as a frame choice. In this paper, we defer making 
a frame-choice until we have derived all the relevant equations, so that we have available 
a set of equations valid for any choice of -u". However, to actually solve the equations, we 
must make a definite choice for -u" (the system of equations is under-determined until such 
a choice is made). Here, it will be convenient to choose u"' to coincide with the velocity of 
the CDM component, since is then geodesic. 

The velocity u°' defines a projection tensor hab which projects into the space perpendic- 
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ular to (the instantaneous rest-space of observers moving with velocity u""): 

hah = gab - UaUb, (2-1) 

where Qab is the spacetime metric. Since hab is a projection tensor it satisfies 

hab = hab), Khcb = hab, K = 3, u'^hab = 0. (2-2) 

We employ the projection tensor to define a spatial covariant derivative (^^V" which 
acting on a tensor T^ ''^d...e returns a tensor which is orthogonal to on every index: 

(3)^a^6...c^ ^ = h^hl . . . h% . . . KWT'-\...^, (2-3) 

where V" denotes the usual covariant derivative. If the velocity field has vanishing 
vorticity (see later) V" reduces to the covariant derivative in the hypersurfaces orthogonal 
to w". 

The covariant derivative of the velocity decomposes as 

WaUb = ^ab + CTab + \Ohab + UaWb, (2-4) 

where Wa = u^VbUa is the acceleration, which satisfies u""Wa — 0, the scalar 9 = V"iia = 3H 
is the volume expansion rate {H is the local Hubble parameter), zuab = ^[aUb] +'W[aUb] is the 
vorticity tensor, which satisfies zUab — T^[ab] and ti^zuab = 0, and aab = ^^^^ (aUb) — dhab /3 is the 
shear tensor which satisfies aab — C(o6)) — ^ '""cTob — 0. The non-trivial integrability 
condition 

= -ZUab<}>, (2-5) 

for any scalar field 0, where an overdot denotes the action of the operator w^Va, follows 
from the Ricci identity. Note in particular that in an evolving universe (0 ^ 0), spatial 
gradients are necessarily non-vanishing in the presence of vorticity. This behaviour, which 
is a consequence of there being no global hypersurfaces which are everywhere orthogonal to 
u"- if the vorticity does not vanish, is central to the discussion of vector perturbations. For 
vanishing vorticity, the 3- Ricci scalar (or intrinsic-curvature scalar) ^^^TZ in the hypersurfaces 
orthogonal to u"- evaluates to 

^^^n = 2Kp-p^ + aaba''', (2-6) 
where p is the total energy density in the u"- frame. 

In an exact FRW universe the vorticity, shear and acceleration vanish identically. In an 
almost-FRW universe, these variables, when suitably normalised to make them dimension- 
less, are regarded as first-order in a smallness parameter e (Maartens et al. 1995). We use 
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the convenient notation, 0{n) to denote that a variable is 0(e"). We assume that products 
of first-order variables can be neglected from any expression in the linearised calculation 
considered here. 

Other first-order variables may be obtained by taking the spatial gradient of scalar 
quantities. Such quantities are gauge-invariant by construction since they vanish identically 
in an exact FRW universe. We shall make use of the comoving fractional spatial gradient of 
the density p^'^ of a species i, 

= ^^'^V„p«, (2-7) 

and the comoving spatial gradient of the expansion 

Z„ = S^^^VaO. (2-8) 

The scalar 5* is a local scale factor satisfying 

S = u''WaS ^HS, (^V"5 = 0(l), (2-9) 

which removes the effects of the expansion from the spatial gradients defined above. The 
vector X^^^ is a manifestly covariant and gauge-invariant characterisation of the density 
inhomogeneity. 

The matter stress-energy tensor Tab decomposes with respect to u"" as 

Tab = PUaUb + ^U^aQb) " phab + T^ab, (2-10) 

where p = Tabu"'u^ is the density of matter (measured by a comoving observer), qa = h^a'^bcU'^ 
is the energy (or heat) fiux and is orthogonal to -u", p = —habT°^l?) is the isotropic pressure, 
and the symmetric traceless tensor i^ab = Kiht'^cd + phab is the anisotropic stress, which is 
also orthogonal to u"'. In an exact FRW universe, isotropy restricts Tab to perfect-fluid form, 
so that in an almost-FRW universe the heat flux and isotropic stress may be treated as first- 
order variables. The final first-order gauge-invariant variables we require derive from the 
Weyl tensor Wabcd, which vanishes in an exact FRW universe due to isotropy. The electric 
and magnetic parts of the Weyl tensor, denoted by £ab ^-nd Bgh respectively, are symmetric 
traceless tensors, orthogonal to u"", which we define by 

£ab = W^W'Wacbd (2-11) 
Bab = —u'Wrjac'^Wefbd: (2-12) 

where rjabcd is the covariant permutation tensor with 770123 = ~\/~g- 
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2.1. Linearised Perturbation Equations for the Total Matter Variables 



Exact equations describing the propagation of the total matter variables (such as the 
total density p), the kinematic variables, and the electric and magnetic parts of the Weyl 
tensor, and the constraints between them, follow from the Ricci identity and the Bianchi 
identity. The Riemann tensor is expressed in terms of Sab, Bab and the Ricci tensor, TZab, 
and the Einstein equation is used to substitute for the Ricci tensor in terms of the matter 
stress-energy tensor. On linearising the equations that result from this procedure (Bruni, 
Dunsby, & Ellis 1992), one obtains five constraint equations: 



Bab + {}^^V'rXJd{a + VVd(„) r76)crU^ 

(P + V)Vab'^U^'^cd + VabcdU^^^'^V^q'^ = 

lK{2^^^VaP + 2eqa + 3^^^V\ab) = 

^'W'zUab + ^'^V'aab-l^'^Vae-nqa = 

^'^V' {Vabcdu''r:u^'') = 

and seven propagation equations: 

Sab + OSab + ^^^V'Bd(^arib)c/u' + |k[3(p + 

+3 ((=^^(,^5) - Ihab^^^^'qc) - 37V,6 - Ql^ab] 








0, 



(2-13) 
(2-14) 
(2-15) 
(2-16) 
(2-17) 



U 



+ -^QOab 



\hab^^^^^Wc) + Eab + ^KTTafc 



Bab + QBab - (^'V^^d(a + V%d(, 

'^ab - ^^^^\aWb] + \Q'^ab 
qa + \Oqa + (p + p)Wa + V^„6 - V^p 
e+\e^- ^^^V'Wa + |/€(p + 3p) 

p + 9{p + p) + ^''^V''qa 









0, 



(2-18) 
(2-19) 
(2-20) 
(2-21) 
(2-22) 
(2-23) 
(2-24) 



where Tab...c = u'^^dTab...c- The constraint equations do not involve time derivatives, and 
so they serve to constrain initial data for the problem. The propagation equations are 
consistent with the constraint equations in the sense that the constraints are preserved in 
time by the propagation equations if they are satisfied initially. The consistency of the exact 
equations follows from their derivation from the exact field equations, and is preserved by 
the linearisation procedure. Including a cosmological constant A in the above equations is 
straightforward; one adds a contribution A/k to the total density p, and subtracts the same 
term from the total pressure p. 



Many of the equations given above have simple Newtonian analogues |( Ellis 1971 )| , and 
thus are simple to interpret physically. The analogues arise because many of the covariantly 
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defined variables have counterparts in (self-gravitating) Newtonian hydrodynamics. An im- 
portant exception is that there is no Newtonian analogue of the magnetic part of the Weyl 
tensor (the electric part is analogous to the tidal part of the Newtonian gravitational poten- 
tial), and no equation analogous to the £ab propagation equation |(Elhs fc Dunsby 1997a) 



These exceptions arise because of the instantaneous interaction in Newtonian gravity, which 
excludes the possibility of gravitational wave solutions to the Newtonian equations. However, 
there is a close analogy between the £ab and Bah propagation equations and the constraint 
equations ( p-14| ) and ( p-15|) , and Maxwell's equations split with respect to an arbitrary 



timelike velocity field (see, for example, Maartens & Bassett (1998)). 

There is some redundancy in the full set of linear equations. (See Maartens (1997) for 
a discussion of the redundancy in the exact non-linear equations for an irrotational dust 
universe.) For example, equation ( ^-13| ), which determines Bab in terms of the vorticity and 



the shear, along with equation ( p-16|) and the integrability condition given as equation ( p-5| ) 
imply equation (|2-14|) . Similarly, equation (|2-19|) follows from equation (|2-13| ) and the prop- 



agation equations for the shear (eq. |P-20|| ) and the vorticity (eq. [P-21|| ). It follows that B, 



■>ab 



may be eliminated from the equations in favour of the vorticity and the shear by making use 
of equation (|2-13|) . This elimination is useful when discussing the propagation of vector and 
tensor modes (see Section |^. 

The usual Friedmann equations describing homogeneous and isotropic cosmological 
models are readily obtained from the full set of covariant equations, since in an exact 
FRW universe the only non-trivial propagation equations are the Raychaudhuri equation 
(eq. ||2-23|| ) and the energy conservation equation (eq. [P-24|| ), which reduce to the Fried- 



mann equation 

H + = -\k{p + ?>p), (2-25) 
and the usual equation for the density evolution 

p = -^H{p + p). (2-26) 

The second Friedmann equation is obtained as a first integral of these two equations: 

H^ + § = l^p, (2-27) 

where 6K/ S'^ is the intrinsic curvature scalar of the surfaces of constant cosmic time. 

The fractional comoving spatial gradient of the density, Xa, and the comoving spatial 
gradient of the are the key variables in the covariant discussion of the 

growth of inhomogeneity in the universe (Ellis & Bruni 1989; Ellis et al. 1989). It is useful 
to have available the propagation equations for these variables. For Xa, we take the spatial 
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gradient of the density evolution equation (eq. ||2-24|| ) and commute the space and time 
derivatives, to obtain 

pXa +{p + p) {Za - SOwa) + Vj3) + Se^^^VaP " OpX^ = 0. (2-28) 
For Za, we take the spatial gradient of the Raychaudhuri equation (eq. | [2-23|| ) which gives 



Za + \eza - s 



- 5(3)Vj=^) V'Wfe = 0, 



(2-29) 



For an ideal fluid (g^ = Ttab = when we choose u"" to be the fluid velocity), with a barotropic 
equation of state p = p{p), the propagation equations for Xa and Za combine with the mo- 
mentum conservation equation (eq. ||2-22|| ) and the integrability condition, given as equa- 
tion ( |2-5| ), to give an inhomogeneous second-order equation for Xa [(Ellis, Bruni, fc Hwang 
1990)|. For a simple equation of state p = (7 — l)p, where 7 is a constant, the second-order 



equation is 



^a+{l- 7) + i(7 - 2)(37 - 2) [19' + ^) Xa 

+ (7 - 1) i^'^v'x^ + 'Sx^ + 27(7 - 1)5^(^) W 



0. 



(2-30) 



From this equation, it is straightforward to recover the usual results for the growth of in- 
homogeneities in an almost-FRW universe [(Ellis fc Bruni 1989 )| . The inhomogeneous term 
describes the coupling between the vorticity and the spatial gradient of the density, which 
arises from the lack of global hypersurfaces orthogonal to u°- in the presence of non-vanishing 
vorticity. In reality, the universe cannot be described by a barotropic perfect fluid. A more 
careful analysis of the individual matter components is required, which we present in the 
next section. 



3. Equations for Individual Matter Components 

In this paper we concentrate on CDM models, so the matter components that we must 
consider are the photons and neutrinos, which are the only relativistic species, and the 
tightly-coupled baryon/electron system and the CDM, which are both non-relativistic over 
the epoch of interest. We consider the description of each of these components separately in 
this section. 
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3.1. Photons 

In rclativistic kinetic theory (see, for example, Misner, Thorne, & Wheeler (1973)), the 
photons are described by a scalar- valued distribution function f^"'\x,p). An observer sees 
f^'^\x,p)(Px(Pp photons at the spacetime point a; in a proper volume (Px, with covariant 
momentum p" in a proper volume (Pp of momentum space. The photon momentum 
decomposes with respect to the velocity it" as 

p" = Eiu" + e^), (3-1) 

where E — p°'Ua is the energy of the photon, as measured by an observer moving with velocity 
It", and Ca is a unit spacelike vector which is orthogonal to u"". 

e"e„ = -1, e-'ua = 0, (3-2) 

which describes the propagation direction of the photon in the instantaneous rest space of the 
observer. With this decomposition of the momentum, we may write the photon distribution 
function in the form f^'^\E, e) when convenient, where the dependence on spacetime position 
X has been left implicit. The stress-energy tensor T^^^ for the photons may then be written 
as 

r}^^ - I dEdn Ef^^ {E, e)paPb, (3-3) 

where the measure dQ denotes an integral over solid angles. The photon energy density p^'^\ 
the heat flux and the anisotropic stress are given by integrals of the three lowest 
moments of the photon distribution function: 

p(T) = J dEdQE^f^^\E,e) (3-4) 
?W = JdEdnE'f^\E,e)ea (3-5) 
ttI? = JdEdnE'f^^\E,e)eae, + lp^^^Kb. (3-6) 

In the absence of scattering, the photon distribution is conserved in phase space. De- 
noting the photon position by x°'{X) and the momentum by p"(A), the path in phase space 
is described by the equations 

p"Vy = 0, (3-8) 

where A is an affine parameter along the null geodesic x^{X). Denoting the Liouville operator 
hy C, we have 

Cf^'\x,p) = i-/(7)(^a(A),p«(A)) = 0, (3-9) 
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in the absence of collisions. Over the epoch of interest here, the photons are not collision- 
less, but instead are interacting with a thermal distribution of electrons and baryons. The 
dominant contribution to the scattering comes from Compton scattering off free electrons, 
which have number density rig in the baryon/electron rest frame. Since the average energy 
of a CMB photon is small compared to the electron mass well after electron-positron annihi- 
lation, we may approximate the Compton scattering by Thomson scattering. Furthermore, 
since the kinetic temperature of the electrons (which equals the radiation temperature prior 
to recombination) is small compared to the electron mass, the electrons are non-relativistic 
and we may ignore the effects of thermal motion of the electrons (in the average rest frame of 
the baryon/electron system) on the scattering. Our final assumption is to ignore polarisation 
of the radiation. Thomson scattering of an unpolarised but anisotropic distribution of radi- 
ation leads to the generation of polarisation, which then affects the temperature anisotropy 
because of the polarisation dependence of the Thomson cross section ax- In this manner, 
polarisation of the CMB is generated through recombination and its neglect leads to errors 



of a few percent |(Hu et al. 1995)| in the predicted temperature anisotropy. We hope to 
develop a covariant version of the radiative transfer equations including polarisation in the 
near future, which should simplify their physical interpretation. 

In the presence of scattering, the photon distribution function evolves according to the 
coUisional Boltzmann equation, 

/:/W(x,p) = c, (3-10) 

where the collision operator for Thomson scattering is 

/i'^)(x,p)-/W(^,P)l, 



C 



(3-11) 



where wi'') is the covariant velocity of the baryon/electron system and f^\x,p) describes 
scattering into the phase space element under consideration: 

3 



IQtt 



/W(a:,p') 



,9 (^V% 



dQ„i(b) , 



(3-12) 



(3-13) 



where e^*-* is the photon direction relative to u^^\ 

Pa = E^'^ {u^:^ + e^:^) , E(^) 

and e'J^^ is the initial direction (relative to m^^-*) of the photon whose initial momentum is 
and final momentum is pa- We write the baryon velocity in the form 

«i^)=7(^)(«« + .i^))' (3-14) 

where vj^'^ is the relative velocity of the baryons, which satisfies u'^vj^'^ = 0, and 7''^'' = 
(1 + g°'^v^''vf'^)~^^'^ . Note that to first-order we have u^^'' = Ua + v^^ since the relative 
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velocities of the individual matter components are first-order in an almost-FRW universe. 
Multiplying the Boltzmann equation by E"^ and integrating over energies, we find 



dEE''Cf^^\E,e 



ne(TT 



-3 



- ne^TT^^^ (l + e^jl"^) J dEE^f^\E, e) 



(3-15) 
(3-16) 



where we have used 

E^') =^^')E{l + e%^^] 

to replace the integral over E by an integral over i?^**^ in the first-term on the right. This 
term can be rewritten as an integral over E''^'^'> = p"^u^^ using the fact that there is no energy 
transfer in Thomson scattering in the rest frame of the scattering electron, so that 



167r 



1 + (g'^^e. 



/(^)(x,p). (3-17) 



(3-18) 



Using the definition of the radiation stress-energy tensor (eq. ( |3-3| )) in the right-hand side, 
we have 

where e^^^ can be expressed as 



7^') (l + e%^J'^ 



(3-19) 



It follows that the energy-integrated Boltzmann equation reduces to 



dEE^Cf'^^\E,e 



rieO-T 



7^'') (l + 



-3 



- UeaT-f^^^ (l + e"t;f ) J dEE^f^\E,e). 



(3-20) 



This equation is exact under the assumption of Thomson scattering and the neglect of 
polarisation. Here, we shall only require the linearised version of equation ( |3-20D ; for a 
covariant discussion of the second-order effects in this equation, see Maartens et al. (1998). 
On linearising equation ( 3-2C ) around an almost-FRW universe, we find 



j dEE^Cf^^\E,e) 



IGtt 



e e 



^(l-4e«.i^))p(^)+.L 

-Hear J dEE^f^\E,e). (3-21) 



This covariant form of the Boltzmann equation was used in Challinor & Lasenby (1998) to 
discuss CMB anisotropics from scalar perturbations on angular scales above the damping 
scale. Note that the equation is fully covariant with all variables observable in the real 
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universe, is valid for arbitrary type of perturbation (scalar, vector and tensor), employs no 
harmonic decomposition and is valid for any background FRW model. 



The numerical solution of the Boltzmann equation (eq. ( p-21| )) is greatly facilitated by 



decomposing the equation into covariantly-defined angular moments. The majority of recent 
calculations (for example, Seljak & Zaldarriaga (1996)) perform an angular decomposition 
of the Boltzmann equation after specifying the perturbation type and performing the ap- 
propriate harmonic expansions. The procedure is straightforward for scalar perturbations in 
a K = universe, where the Fourier mode of the perturbation in the distribution function 
may be assumed to be axisymmetric about the wavevector k (this assumption is consis- 
tent with the evolution implied by the Boltzmann equation), allowing an angular expansion 
in Legendre polynomials alone. However, for tensor perturbations the situation is not so 
straightforward (see, for example, Kosowsky (1996)), since the Boltzmann equation does 
not then support axisymmetric modes. Instead, the necessary azimuthal dependence of the 
Fourier components of the perturbation in the distribution function, which is different for the 
two polarisations of the tensor modes, must be put in by hand, prior to a Legendre expansion 
in the polar angle. This procedure may be eliminated by performing a covariant angular 
expansion of f^'^\x,p) prior to specifying the perturbation type or background FRW model. 
The covariant (tensor) moment equations that result may then be solved for any type of per- 
turbation (and any background curvature K) by expanding in covariant tensors derived from 
the appropriate harmonic functions (see Section ^ for the case of scalar perturbations and 
Section ^for tensor perturbations). This procedure automatically takes care of the required 
angular dependencies of the harmonic components of the distribution function, allowing a 
streamlined and unified treatment of all perturbation types in background FRW models with 
arbitrary spatial curvature. 

The covariant angular expansion of the photon distribution function takes the form 
(Ellis, Matravers, & Treciokas 1983, Thorne 1981) 

oo 

f^^\E,e)=J2F^[l,^e'^^e^^...e^', (3-22) 
/=o 

where the tensors have an implicit dependence on spacetime position x and energy 

E, and are totally symmetric, traceless and orthogonal to m": 

^ p(0 „aia2 = n yO-i f(0 = n ('3-23') 

-^ai...ai -^(ai...ai)' aia2...ai '-'5 " ai...ai I'-' ^'-'J 

Employing the expansion given in equation ( |3-22|) , the action of the Liouville operator on 
f^^\E,e) reduces to 
00 

/:/W(E, e) = E [dEFlll^^e^^d^E + /V,Fi;)..,,e"^ + /F^^^/V^e"^] e'^^ . . e'^K (3-24) 
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Using the geodesic equation, we find that 

hahP^^ce^ = -E (abce^e^ea + aahe^ + w'^CaCb + Wa- tUabe^) , (3-25) 



which is first-order. In an exact FRW universe, isotropy restricts Fj^^ ai — ^ f*^^ / > 0, so 
that in an almost-FRW universe 

^al.ai = 0(1) for / not equal to zero. It follows that the 
last term in equation (|3-24|) makes only a second-order contribution and may be dropped in 
the linear calculation considered here. 

Inserting the expansion given in equation ( |3-22D into the Boltzmann equation (eq. [ ^-21[ |) 
and performing a covariant angular expansion of the resulting equation gives a set of moment 
equations which are equivalent to the original Boltzmann equation. The linearised calculation 
is straightforward, although a little care is needed for the first three moments since F^^^ is a 
zero-order quantity. (The exact expansion of the left-hand side of the Boltzmann equation, 
equation ( P-24| ), is given in Ellis et al. (1983) and Thorne (1981).) For / = 0, 1 and 2, we 
find 

^i^^ + + VV« + |pW^, - |(=^)V,p(-) = n^ar (fp^^^t;!^) - gi^^) (3-27) 

= -^^n,aTni^\ (3-28) 

and, for / > 3, 

•Jai...a, ^ 3^'Jai...ai ^ ^ -'bai...a, (21+1) \ ^ (^i "'aa-.a,) 

The tensors J^'-* which are traceless, totally symmetric and orthogonal to u"", are derived 
from the Fj^^ by integrating over energy: 

^a....a. - (2/+l)(2/)!yo ^^^^^ 

The constant factor is chosen to simplify algebraic factors in the moment equations. Using 
equations ( P-4| - |3-6| ), the lowest three moments relate simply to the energy density, heat flux 
and anisotropic stress: 

pW = j(^ gW = Ji^), vr(? = ji?. (3-31) 
It is straightforward to show that the tensor 
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which appears in equation ( P-29|) is traceless, symmetric, and orthogonal to as required. 

It will be observed that for / > 3, the moment equations link the / — 1, / and / + 1 
angular moments of the (integrated) distribution function, while the / = 2 equation also 
involves the density p^'^'> which is the / = moment. The exact moment equations that 
arise from expanding the Liouville equation in covariant harmonics also couple the / + 2 
and / — 2 moments to J^'-* (Ellis et al. 1983), but these terms are second-order for / > 3 
and so do not appear in the linearised equations presented here. In the exact expansion of 
the Liouville equation, the coefficient of the / + 2 angular moment in the exact propagation 
equation for J^'^^ is the shear aab-, which leads to the result that the angular expansion of 
the distribution function for non-interacting radiation can only truncate {J^^ a, = 
I greater than some L) if the shear vanishes |( Ellis 1996)| . This exact result, which is lost in 



linearised theory which permits truncated distribution functions with non-vanishing shear, 
is an example of a linearisation instability (see Ellis & Dunsby (1997b) for more examples). 
However, this is not problematic for the linearised calculation of CMB anisotropies since it 
is never claimed that the higher-order moments of the photon distribution vanish exactly. 
Instead, the series is truncated (with suitable care to avoid reflection of power back down 
the series) for numerical convenience. The truncation is performed with L large enough so 
that there is no significant effect on the J'^^ for the range of / of interest. 

Finally, by taking the spatial gradient of equation (|3-26|) , and commuting the space 



and time derivatives, we find the propagation equation for the comoving fractional spatial 
gradient of the photon density, 

+ |Z„ + -^(^)v^(^)v'qi^^ - ISOwa = 0, (3-33) 

where Za is the comoving spatial gradient of the volume expansion. 



3.2. Neutrinos 

We consider only massless neutrinos, and these are non-interacting over the epoch of 
interest. It follows that their distribution function f^'^\x,p) satisfies the Liouville equation 
Cf^^\x^p) = 0. Expanding the neutrino distribution function in covariant angular harmon- 
ics, we arrive at the moment equations for the tensors G^J:^ , which are defined in the same 
manner as the J^^^ but with the photon distribution function replaced by the neutrino 
distribution. These moment equations are the same as the photon equations, but with the 
scattering terms omitted: 
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pH + l^pH + (3)v«gH = (3-34) 

^i'^^ + l^^i'^^ + ^'^V^vri^ + |p('^)-|(3)vvH = (3-35) 

^i? + l^vri^ + V^GSi - I ((^)V(.g(^) - IhJ'^V^gi^^) - ^p^^V., = 0, (3-36) 

and for / > 3, 

^ai...ai^ 3'^'^ai...ai+ ^^bai...ai (2l+l)\ ^ [ai^ a2...ai) 

- M^^^V^Gg-L,_.^.-..)) = 0. (3-37) 

The propagation equation for the comoving fractional spatial gradient of the neutrino density, 
follows from equation ( P-34|) : 

^i"^ + l^a + ^T^'V'^^^V^gi^) - 1^^^. = 0. (3-38) 



3.3. Baryons 

Over the epoch of interest here, the electrons and baryons are non-relativistic, and 
may be approximated by a tightly-coupled ideal fluid (the coupling arising from Coulomb 
scattering). The energy density of the fluid is p^''\ which includes contributions from both 
the baryonic species and the electrons, the fluid pressure is p^^\ and the velocity of the fluid 
is u^^"* = Ua + v^Ji'^ to first-order, where the 0(1) relative velocity v^^ satisfies u"'vj^^ = 0(2). 

The linearised baryon stress-energy tensor evaluates to 

T^^ = p^%aut - P^'^K, + 2(p(^) + p^'^)u^avt^, (3-39) 

which shows that there is a heat flux (p(^) +p(''))t;(^) due to the baryon motion relative to the 
M° frame. The equations of motion for p*^''^ and v^"* follow from the conservation of baryon 
plus photon stress-energy (the baryons and photons interact through non-gravitational effects 
only with themselves): 

v"r„?) + v^ri^) = 0. (3-40) 

Using the / = and / = 1 moment equations for the photon distribution, we find the 
propagation equation for the baryon energy density: 



(3-41) 
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and a propagation equation for v^^ : 



+ n^ar {Ip^M'^ - gi^^) = 0, (3-42) 



which must be supplemented by an equation of state hnking p^^^ and p^^K The final term 
in equation ( P-42| ) describes the exchange of momentum between the radiation and the 
baryon/electron fluid as a result of Thomson scattering. There is no such term in equa- 
tion ( p-41| ) since both the radiation drag force and the baryon velocity relative to the m° 
frame are first-order, which give only a second-order rate of energy transfer in the frame. 
Energy transfer due to thermal motion of the electrons in the baryon rest frame has negli- 
gible effect on p^^^ since the electrons are non-relativistic; ksT^'^'^ <C me, where T*^^^ is the 
baryon kinetic temperature (assumed equal to the electron kinetic temperature), and rUf. is 
the electron mass. 



Taking the spatial gradient of equation (|3-41D gives the propagation equation for 



the fractional comoving spatial gradient of the baryon energy density: 

+ Se^^WaP^'''^ - Op^^'^X^^^ = 0. (3-43) 

We have retained all terms involving the baryon pressure p^'^"' in the equations of this section. 
In practice, over epochs where the baryons are non-relativistic {p^^^ -C p^^^), the only pressure 
term that need be retained is the term (3)V"p*^^^ which appears in equation ( P-42| ) . This term 
appears as a small correction to the total sound speed in the tightly-coupled baryon/photon 
plasma, and is potentially significant during the acoustic oscillations in the plasma. 



3.4. Cold Dark Matter 

We will only consider cold dark matter (CDM) here, which may be described as a 
pressureless ideal fluid. Hot dark matter (HDM), which requires a phase space description, 
is considered in Ma & Bertschinger (1995), where the calculations for scalar perturbations 
are performed in the synchronous and conformal Newtonian gauges. The CDM has energy 
density p^'^^ in its rest frame, which has velocity u*^^ = Ua + with the first-order relative 
velocity v^f^ satisfying u"'v^^ = 0(2). The linearised CDM stress-energy tensor evaluates to 

r„^) = pWw„w, + 2pWw(,t;Jj\ (3-44) 
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which is conserved since the CDM interacts with other species only through gravity. The 
conservation of T^^^ gives the propagation equations for p'^'^^ and v'^'>: 

v^^ + \dv^^ + Wa = 0. (3-46) 

Since the CDM moves on geodesies, the velocity m^^^ provides a convenient frame choice. With 
this choice, the acceleration Wa vanishes. We use the CDM frame to define the fundamental 
velocity in Section ^, where we discuss scalar perturbations in the CDM model. For 
the moment, however, we continue to leave the choice of frame unspecified for generality. 
The final equation that we require is the propagation equation for the fractional comoving 
spatial gradient of the density, This follows from equation ( P-45|) on taking the spatial 

gradient: 

A'i'^) + Za + S^'^^Vj^^V'vi'^ - S9wa = 0. (3-47) 

The equations for the matter components that we have described in this section combine 
with the covariant equations of Section ^ to give a complete description of the evolution 
of inhomogeneity and anisotropy in a fully covariant and gauge-invariant manner. The 
equations given in Section |^ make use of the total energy density and pressure, heat flux and 
anisotropic stress. These quantities are related to the individual matter components in the 
CDM model by 

•,(7) I ni^) ^ nib) , (c) 



p = p(^) + pW + pW + pC^) (3_48) 
p = ip(7) + ipH+pW (3.49) 

qa = q^:^ + q^:^ + {p^'^ + P^'^) v^:^ + p^'^hi^^ (3-50) 

= vri^+vri^. (3-51) 

The equations given here are both covariant and gauge-invariant. Employing gauge- 
invariant variables ensures that the problem of gauge-mode solutions does not arise, and 
that all quantities are independent of the choice of map between the real universe and a 
background FRW model. We have only considered the linearised equations here, but the 
linearisation procedure is not fundamental to the covariant and gauge-invariant approach. 
It is straightforward to extend the treatment to include non-linear effects (Maartens et al. 
1998), which should provide a systematic footing for the discussion of second-order effects 
on the CMB. Indeed, the simplicity with which exact, non-linear equations can be written 
down and manipulated is a significant virtue of the covariant approach. Unlike in Bardeen's 
gauge-invariant approach |(Bardeen 1980)| , the definition of the variables employed here does 



not require that the perturbations be in the linear regime, and furthermore, the variables 
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do not depend on the non-local decomposition of the perturbations into scalar, vector and 
tensor type and the associated harmonic analysis. The covariant approach describes scalar, 
vector and tensor modes in a unified manner, although decomposing the linear perturbations 
is useful to aid solution of the linearised equations late on in the calculation. A further 
advantage of the covariant and gauge-invariant approach over that introduced by Bardeen, is 
that only covariantly defined variables are employed, which are simple to interpret physically. 
In contrast, the Bardeen variables are constructed by taking linear combinations of (gauge- 
dependent) metric and matter perturbations in such a way that the resulting variable is 
gauge-invariant (for small gauge-transformations which preserve the scalar, vector or tensor 
structure of the metric perturbation). These variables have simple physical interpretations 
only for certain specific gauge choices. Finally, note that we have not yet had to specify 
whether the background FRW model is open, flat or closed. However, we have made the 
implicit assumption that the universe is almost FRW when specifying the zero and first-order 
variables in the linearisation procedure. 



4. The CMB Temperature Anisotropy 

The energy-integrated moments J^^ of the photon distribution function provide a 
fully covariant description of the CMB temperature anisotropy. In the u"" frame, denote the 
average bolometric temperature on the sky at the point x by To(x), so that 

oc i-| dEdnE'f(^\x,p) = i-jW, (4-1) 

which is just the Stefan-Boltzmann law. We use the fractional temperature variation 5r(e) 
from the full-sky average Tq to characterise the temperature perturbation along the spatial 
direction e in a gauge-invariant and covariant manner (Maartens et al. 1995; Dunsby 1997). 
It follows that 

(l + 5He))' = / dEE'f^^\x,p), (4-2) 

so that to first-order 

The right-hand side of equation ([4-3|) is the covariant angular expansion of the temperature 
anisotropy. The tensors Ja}...ai ^^us provide a natural covariant description of the CMB 
anisotropy. They may be related to the more familiar aim components in the spherical 
harmonic expansion of 5r(e) by introducing an orthogonal triad in the instantaneous rest 
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space at x, so that e"" = (sin^ cos0, sin^ sin0, cos^). Then the two expansions are related by 

E ^ a^MO, ^) = ^^(^J^y (ff.^ ' 4?...a.e" ■ ■ ■ e"'. (4-4) 



1 5Z ^4^(,)N2(_2)^f/!)2^'^''''---^"'''4^...a,4?..V (4-5) 



Squaring this expression and integrating over sohd angles, we find the important rotational 
invariant 

1 ^ , ,2 _ 47r (2/)! 
2/ + l,iir'"^' - (4p(7))M-2)'(/!) 

The quantity on the left is a quadratic estimator for the CMB power spectrum (the Ci), which 
we see is related to the covariant tensors J^'^ in a very simple manner. Further properties 
of the covariant and gauge-invariant description of CMB temperature anisotropies are given 
by Gebbie & Ellis (1998). 



5. Scalar Perturbations 

Up to this point, we have treated the scalar, vector and tensor modes of linear theory 
in a unified manner. However, to obtain solutions to the covariant equations it proves 
useful to consider scalar, vector and tensor modes separately. In this section we consider 
scalar modes; tensor modes are discussed briefiy in Section |^. (Vector modes decay in an 
expanding universe in the absence of defects, and so are not likely to have a significant 
effect on the CMB in infiationary models.) In the covariant approach to perturbations in 
cosmology, we characterise scalar perturbations by demanding that the magnetic part of the 
Weyl tensor and the vorticity be at most second-order. Setting Bab = 0(2) ensures that 
gravitational waves are excluded to first-order, and demanding that vjab = 0(2) ensures 
that the density gradients seen by an observer in the frame arise from clumping of the 
density, *^^^V^p = 0(1), and not from kinematic effects due to vorticity (the absence of 
flow-orthogonal hypersurfaces), which give "^^^V^p = 0(2) in an almost-FRW universe. Note 
that we do not classify scalar perturbations as having Bab = ^ab = (to all orders), which 
is only a highly restricted subset of the full set of scalar solutions. For example, in an 
(exactly) irrotational dust-filled universe (a "silent" universe), it can be shown from the 
exact non-linear equations that demanding Bab = forces the solution into a very small 



class, which probably all have high symmetry [(Ellis 1996 )| , and so cannot represent a very 
general perturbation. This arises because requiring that Bab = be preserved along the fiow 
lines introduces a series of complex constraints which reduce greatly the size of the solution 
set. However, requiring only that Bab and Wab be at most second-order gives a much larger 
class of solutions because only two new constraints are introduced, and these are necessarily 
preserved by the propagation equations. The solutions with Bab and zuab vanishing exactly 
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comprise a very small subset of the larger class of exact solutions which we classify as scalar 
perturbations. 

On setting Bab = and zuab = 0, (equality to zero in the linearised theory should be 
taken to imply that the quantity is at most second-order) we see from equation (|2-13|) that 

VV,(,r/,),/M^ = (=^)V[,(=^) VV,]e = 0, (5-1) 

where the antisymmetrisation is on the indices a and b in the right-hand equation. This 
is a necessary condition for aab to be constructed from a scalar potential. It follows from 
equations (p-16| ) and (|2-21D that 

^'^VlaQb] = 0, ^^^VlaWb] = 0, (5-2) 

so that the heat flux and acceleration may be written as spatial gradients of scalar fields 
(making use of the integrability condition given as equation ( |2-5| )). Consistency of ^^^V[aqb] = 
with equation ( ^-22| ) for then requires that 

^')V[„(3)V%,], = ^ ^'^V[a^^^V''£b]c = 0, (5-3) 

with the implication following from equation (|2-15|) . It follows that all vector variables, such 
as Qa and ^^■'V^'^ab, may be derived from scalar potentials. The new constraint, given as 
equation ( |5-1| ), is only consistent with the propagation equations if Sab and Hab satisfy 

^'^V^£:,(,r/6)ee"w^ = ~l^'^'^V'MaVb)ce''u'- (5-4) 



In the absence of anisotropic stress, we see that the left-hand side of equation (5^) is 
constrained to be zero, which is consistent with the propagation equation for Sab, given as 
equation ( |2-18| ), with Hab = 0. If the anisotropic stress does not vanish, we include the 
constraint 

'-^^V''SdiaVb)ce''u' = ('V%,(ar/,),/M^ = 0, (5-5) 

in the definition of a scalar mode, which is easily shown to be consistent with the propagation 
equation for Sab- Requiring consistency of equation ( |5-5| ) with the propagation equation for 
TTab implied by the photon and neutrino Boltzmann hierarchy yields a series of constraints 
on the moments Ja}...ai and G^aj a^, which are necessary conditions for them to be derived 
from scalar potentials. 

The new constraint equations that we have introduced, by restricting the solution to be 
a scalar mode, may be satisfied by constructing the covariant and gauge-invariant variables 
from tensors derived from scalar potentials by taking appropriate spatial covariant derivatives 
of the scalar functions. It proves convenient to separate the temporal and spatial aspects of 
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the problem by expanding the scalar potentials in the eigenf unctions Q^^^ of the generalised 
Helmholtz equation (Hawking 1966; Ellis et al. 1989) 

(3)v2qW ^ {3)v"(3)v,Q('=) = ^Q^^\ (5-6) 

which are constructed to satisfy 

Q^'^) = 0. (5-7) 

These equations are correct to zero-order only; in such equations, equality should be under- 
stood to this order only. In general, we cannot impose equation ( ^-6|) and Q^^-* = exactly, 
since the constraint equation is inconsistent with the evolution implied by Q'^^'^ = at first- 
order. The allowed values of the eigenvalues k'^ / S"^ are determined by the scalar curvature 
of the background model (since Q'^^'^ are only needed to zero-order). In a flat model, = 0, 
is a comoving continuous wavenumber > 0. In closed models, K > 0, k takes only discrete 
values with k'^ = 7(7 + 2)7^ where 7 is a non-zero, positive integer. In open models, K < 0, k 
again takes continuous values, but with the restriction k"^ > \K\. More details may be found 
in Harrison (1967). The eigenfunctions Q^'^^ are labelled by the lumped index k. This index, 
which determines the eigenvalue k'^/S'^, should be understood to distinguish implicitly the 
distinct degenerate eigenfunctions which all have the same eigenvalue k'^/S'^. This multiple 
use of the symbol k should not cause any confusion, since the lumped index will always 
appear as a superscript or subscript. A function of the eigenvalue k will be denoted with the 
eigenvalue as an argument, for example A{k), to distinguish it from the quantity which 
depends on the mode label k and not just the eigenvalue. From the Q*-*'-' we form a vector 

Ql^^ ^ l'^'^VaQ^'\ (5-8) 
which is orthogonal to u"' and is parallel transported at zero-order along the flow lines: 

W'Qi'^ = 0, Qi'^ = 0. (5-9) 

We define totally symmetric tensors of rank /, Q^^-* by the recursion formula (for / > 1) 
(see also Gebbie & Ellis (1998)) 

r,{k) _ 5 A3)y r,{k) _ ±i(3)yb-o(fc) \ (r,TQ\ 

Vai...a, - fc V ^K^a2...a0 21-1 ^ ^b{ai...ai_2'^ai_iai) ) ■ {^^^) 

These tensors satisfy the zero-order properties 

«"^Qil...a, = 0, /^"^'^^Qil...., = 0, = 0, (5-11) 

which are readily proved by induction. 

The scalar functions Q^^^ are the covariant generalisations of the scalar eigenfunctions of 
the Laplace-Beltrami operator on the homogeneous spatial sections of the background FRW 
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model, which are usually employed in the harmonic decomposition of perturbed quantities 
(see, for example, Bardeen (1980)). In the covariant approach, attention is focused on 
a velocity field ■u'*, rather than a spatial slicing of spacetime, so it is natural to employ 
harmonic functions defined by equation (5^). Some of the differential properties of the 
derived tensors Q^ai...ai (fo^ ^ ^ 2) are given in the appendix to Bruni et al. (1992). We add 
two more results to this list which will be useful later: 



(3)v7ai/n{fc) 

* ^ aia2...ai 



(3)y2g(A0 



I 



21-1 S 

hi 

52 



1) 



fc2 



a2...ai 



)(fc) 



(5-12) 
(5-13) 



which may be derived from the recursion relation given as equation ( |5-10|) and the definition 



in equation ( |5-6|) . Further properties of the scalar harmonics are given by Gebbie & Ellis 
(1998). 

To relate the timelike integration employed in the Boltzmann multipole approach. 



adopted here, to the lightlike integration employed in line of sight methods (Seljak fc Zaldar- 



|riaga 1996)| , which we consider in Section ^.2| , it is convenient to note the following zero-order 
results for the variation of the Qi^^ ai along the line of sight through some point R. Let x'^(A) 
be a null geodesic with tangent vector parallel to Ua + Cq, and A satisfying {ua + Ca) VA = 1. 
Define a positive parameter y{\) along the past null geodesic by dy/d\ = —k/S with y = 
at the point of observation R. Then the evolution of the quantities Q^^ „^e"^ . . . e"' is given 
to zero-order by the hierarchy 
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which follows from the recursion relation given as equation (|5-10|) . We shall only consider the 
solution to this hierarchy in a = universe, in which case we find the following variation 
of Q^*^-* along the line of sight: 



(Q 



^ ^-2)'(/!)2 ^^^y^^^^<^i- 



ai 



(5-15) 



where the ji{y) are spherical Bessel functions. Equation (|5-15[ ) expresses Q'^'^^ a parameter 



distance y down the line of sight in terms of the Q^^-* ^^e"^ . . . e"' at the point of observation. 



R. If required, the value of Q^^^ ^^e"^ . . . e'*' down the lightcone can be found from the solution 
for Q(^) (eq. ( |5-15|) ) and the hierarchy (eq. (|5-14|) ). 

The additional constraints introduced by the conditions for a scalar mode are satisfied 



- 26 - 



identically if we construct the gauge-invariant variables in the following manner: 
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Wa 




(5-23) 
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where i labels the particle species (and we omit the label when referring to total fluid 
variables). The symbolic summation in these expressions is a sum over eigenfunctions of 
equation ( |5-6| ). For closed background models, the sum is discrete, but in the fiat and open 
cases the summation should be understood as an integral over the continuous label fc, which 
distinguishes distinct eigenfunctions. The scalar expansion coefficients, such as \ are 
themselves first-order gauge-invariant variables, which satisfy 

(3)^,^« = 0(2). (5-24) 

They are labelled by the lumped index k. Finally, we assume that the higher-order angular 
moments of the photon and neutrino distribution functions may also be expanded in the 
Q^^^ harmonics. By considering the zero-order form of the scalar harmonics Q^'^\ and 
derived tensors, it is straightforward to show that this condition is equivalent to the usual 
assumption that the Fourier components of the distribution functions are axisymmetric about 
the wavevector k (see, for example, Seljak (1996)). With this condition, we have 

Jt.a,=p'^'i:4'Q'S...a, G^U=p^^^Y.Gi:^Qil,^, (5-25) 

k k 

for photons and neutrinos respectively. 
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5.1. The Scalar Equations 



It is now a simple matter to substitute the harmonic expansions of the covariant variables 
into the constraint and propagation equations given in Sections |^ and ^ to obtain equations 
for the scalar expansion coefficients that describe scalar perturbations in a covariant and 
gauge-invariant manner. To simplify matters, we assume that the variations in baryon 
pressure p'-^-' due to entropy variations are negligible compared to those arising from variations 
in p'^'^\ so that we may write 



(5-26) 

where is the adiabatic sound speed in the baryon/electron fluid (this is different to the 
total sound speed in the tightly-coupled baryon/photon fluid). 

With this assumption, we obtain the following equations for scalar perturbations: for 
the spatial gradients of the densities, we find 



(7) 



(7) 



X, 
X, 



■Ac) 



ib) 



S \3^' 

-I + ) + lew, 

(1 + 



9wk 



+ {^- cl) ex, 



k ! 



(5-27) 
(5-28) 
(5-29) 
(5-30) 



for the spatial gradient of the expansion, we find 



2 (pW^-) + p(^)^i^)) + p(^)A'f + (1 + ^cDp^^xi 
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|kp7 - ^ 



(5-31) 



where 7 is defined in terms of the total pressure p and density p by p = (7 — l)p (note that 
we do not assume that 7 is constant). The heat fiuxes satisfy 

= Hear (I? 
= 0, 



ql:^ + lU2nr-Xr+^w,) 
and for the baryon and CDM peculiar velocities 

+ l^^k + -Q^k 



(5-32) 
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The propagation equations for the anisotropic stresses are 
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and the remaining moment equations, for / > 3, are 
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The propagation equations for £ab and aah become 

(ff + \0^k) + \^^^p{,l(7k + qk) + |«;p^(37 - l)vrfc - i/tpTrfc 

+ l^dfc) + (1)^ {^k - Wk) + \^^p'Kk 
Finally, the remaining constraint equations become 



(0 




0. 



2(|f(l-f)$. 



m 



Xk 4- 
Zk 



1 3ft: \ ^ 

1 - TTfc 



1 



n- 



- K,pOqk 
+ i^PQk 




0. 



(5-38) 
(5-39) 

(5-40) 
(5-41) 

(5-42) 
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The variables Xk-, Qk and vr^ refer to the total matter, and are given in terms of the component 
variables by 



P^k 

pqk 



(5-44) 
(5-45) 
(5-46) 



These equations give a complete description of the evolution of inhomogeneity and 
anisotropy from scalar perturbations in an almost-FRW universe with any spatial curvature. 
The system closes up once a choice for the velocity is made, and it is straightforward 
to check that the constraint equations are consistent with the propagation equations. The 



equations for J^^^ and G^^' for Z > 3 are equivalent to those usually found in the literature 
(see, for example. Ma & Bertschinger (1995) and set K = 0). This is because the moments of 
the perturbed distribution function, used in such gauge-dependent calculations, are gauge- 
invariant for / > 1. (The / = 1 moment does depend on the choice of coordinates in the 
real universe, but is independent of the mapping onto the background model, since the 
background distribution function has no angular dependence.) Gauge- invariant versions of 
the usual synchronous-gauge equations |(Ma fc Bertschinger 1995)] are obtained by taking m° 
to coincide with the CDM velocity, so that Wa and v^^ vanish. 



^(0 



5.2. The Integral Solution 



The integral solution to the Boltzmann multipole equations is central to the line of 
sight integration method for the calculation CMB anisotropics. This method, which has been 
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implemented in the CMBFAST code of Seljak & Zaldarriaga (1996), provides a very fast route 
to calculating the CMB power spectrum. Although we do not make use of this method for 
the numerical calculations presented in this paper, it may be useful to some readers to have 
available the integral solution to the covariant and gauge-invariant Boltzmann hierarchy, not 
least because it provides the link between the lightlike integrations along the observational 
null cone, employed in line of sight methods (which includes the original calculation by Sachs 
& Wolfe (1967)), and the timelike integrations along the flow lines of the velocity field u", 
employed in the Boltzmann multipole approach. For simplicity, we restrict attention to 
K = almost-FRW universes. 

The (linearised) integral solution to the hierarchy of Boltzmann multipole equations 
given in the previous subsection is, for / > 1, 



R 



tR 



(5-47) 



where ( )r denotes the quantity evaluated at R, and the integral is taken along the flow line 
of u"" through the point R. Here, t is proper time along the flow line (with t = tji a.t R), 
y = y{tR,t') is k times the conformal time difference along the flow line between t' and 1^: 
y = J^,^ k/Sdt, and f = f{tR,t') is an "optical depth" along the flow line: f = J^,'^ rieaT dt. 
In deriving equation (|5-47|) we have neglected any contribution from initial conditions, since 
these are exponentially suppressed by a factor exp[—f (tji, 0)]. It is straightforward to verify 
by differentiating with respect to that equation ( |5-47D is the solution to the Boltzmann 
hierarchy for scalar perturbations in a i^' = almost-FRW universe. Verification for the 



/ = 1 moment requires the following formal solution for Xj: 



(7). 
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tR 
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iMy) + j^kiy) 
3o{y)}dt', 



po{y) 

(5-48) 



where again we have neglected the exponentially suppressed contribution from the initial 
conditions. In numerical applications, it is convenient to manipulate equation ( 5-47| ) further 
by integrating by parts |( Seljak fc Zaldarriaga 1996) . 



The integral in equation ( |5-471 ) is taken along the flow line of u"' through R. However, 
in the linearised calculation considered here, the integral can be performed along (any) null 
geodesic through R also. To see this, regard y and f as the restrictions to the flow line 
of (zero-order) fields in the past lightcone of R, with ^^^Vajj = 0(1) and similarly for f. 
Replacing the measure dt' by Uadx"" in the integral in equation ( |5-47| ) , and noting that 
V[aMfe] = 0(1) and, for example, *-^^VaCrA; = 0(2), it follows that the line integral of times 
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the integrand in equation ( 5-47| ) is path-independent. The flow hne and null geodesic through 
R can be joined at early times by a spacelike curve with Uadx"" = 0(1), which therefore makes 
only a second-order contribution to the integral around the closed loop. To zero-order, the 
restriction of the fields y and f to the null geodesic through R define quantities y = y{XR, A') 
and r = r{XR,X') on the null curve, where y = Jy^k/SdX and r = J^'^ rigaT dX, so that r 
is the optical depth along the line of sight. (The parameter A along the null geodesic was 
defined in Section |^.) Using the integral along the line of sight, we can now reassemble the 
gauge-invariant temperature perturbations from the mean, Srie), at R using equations (13) 
and ( |5-25|) . Recalling equation ( |5-15|) for the variation of the quantities Q^ai...ai^°'^ . . . e""' down 
the line of sight, the temperature anisotropy from scalar perturbations in an almost-FRW 
universe reduces to 



(Srie)) 
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^vr^Qiye^e^ - + \4'^Q^'^ dX'. (5-49) 



Equation ( |5-49|) , first given in this covariant form in Challinor & Lasenby (1998), is valid for 
any value of the spatial curvature, even though the derivation given here considers the K = 
case only. Equation ( |5-49| ) is most easily obtained by a direct integration of the covariant 
Boltzmann equation for the temperature anisotropy Sxie) along the line of sight [((Jhallinor 
fc Lasenby 1998)| . However, the route followed here makes clear the link between the lightlike 
and timelike integrations, employed in the line of sight and Boltzmann multipole methods 
respectively. 



5.3. Initial Conditions on Super-Horizon Scales 

In this subsection, we analytically extract the solution of the scalar perturbation equa- 
tions in the radiation dominated era. We shall only consider modes with \K\/k'^ <C 1 so that 
we may ignore terms involving K in the scalar equations. Associated with each mode there 
is a characteristic length scale, S/k. The condition \K\/k'^ <^ 1 is equivalent to requiring 
that this length scale be small compared to the curvature radius of the universe. For such 
modes, k is effectively a comoving wavenumber. We shall also require that the mode be well 
outside the horizon scale l/i?, so that we consider only those modes satisfying 



(5-50) 
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where Tik = SH/k is the ratio of the characteristic length scale to the horizon scale, and 
H'^S'^/\K\ is the (squared) ratio of the curvature radius to the horizon scale. If the universe 
may be approximated by a i^" = universe to zero-order, equation ( |5-5(JD reduces to Tik ^ 1- 
The approximate analytic solution may be used to provide initial conditions for a numerical 
integration of the scalar equations (see Section P). 

At this point, it is convenient to make a choice of frame. In the CDM model, the rest 
frame of the CDM defines a geodesic frame, which provides a convenient choice for m", since 
the acceleration then vanishes identically. We assume that this frame choice has been made 
in the rest of this paper. 

Well before decoupling, the baryons and photons are tightly-coupled because of the high 
opacity to Thomson scattering. This scattering damps the photon moments for I > 2, but 
a dipole (/ = 1) moment can survive if the baryon velocity does not coincide with the CDM 
velocity. To a good approximation, we may ignore the J^^ for I > 2, and set v^^^ = Sq^"^ /A 
so that the radiation is isotropic in the rest frame of the baryons. This is the lowest-order 
term in the tight-coupling approximation (see Section |6.2| ). Similarly, we expect that the 
higher-order neutrino moments will also be small in the early universe, since the neutrinos 
were in thermal equilibrium prior to their decoupling. Furthermore, the baryon and CDM 
densities, p^''^ and p^^\ are negligible compared to the radiation and neutrino densities, p^"'^ 
and p(^), in the radiation dominated era. 

A useful first approximation to the full set of scalar equations is obtained by setting the 

.(0 



neutrino moments, G). , to zero for I > 2. It is convenient to take the dependent variable to 
of tl 

equations of the previous section reduce to the following set: 



be X = ^ instead of the proper time t along the flow lines, so that the scalar propagation 
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(5-51) 

(5-52) 
(5-53) 
(5-54) 
(5-55) 
(5-56) 
(5-57) 



where a prime denotes differentiation with respect to a;, and we have used the zero-order 
Friedmann equation in the form = Kp/3 since the curvature term may be neglected by 
equation (|5-50D. We have followed Ma & Bertschinger (1995) by introducing the dimension- 
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less quantity R defined by 



R = 



P 



(5-58) 



After neutrino decoupling, i? is a constant which depends only on the number of neutrino 
species. The remaining equations that we require are the two scalar constraints which reduce 
to 



2x^$A. - 3a; 
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(5-59) 
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This set of equations gives a closed equation for 

3x^1 + 12$'fc + x^k = 0. (5-61) 

This equation should be compared to the fourth-order equation for the metric perturba- 
tion variable in the synchronous gauge (see, for example. Ma & Bershinger (1995)). The 
fourth-order equation admits four linearly-independent solutions, but two of the solutions 
are gauge modes that arise from mapping an exact FRW universe to itself. The gauge- 
invariant approach adopted here ensures that such gauge modes do not arise. This is evident 
from equation (^-6 IJ ) which is only a second-order equation. The two linearly-independent 
solutions of this equation both describe physical perturbations in the Weyl tensor, which 
vanishes for an exact FRW universe. It is now straightforward to find the general solution of 
equations ( |5-51[ )-( p-6CI|) . There are two solutions with non- vanishing Weyl tensor ($fc 7^ 0), 



(7) 



which we write as 

= -3y'^[iCy-D)cosy-iC + Dy)smy] (5-62) 

Zk = 3V3y~^[2{C + Dy)cosy + 2{Cy-D)smy-C{2 + y^)] (5-63) 

sVSy-^ [D cosy + C siny - Cy] (5-64) 

-4:VSy-^[C cosy - D siny - C] (5-65) 

~^y-^[{2RC + Dy) cosy + {Cy - 2RD) siny - 2RC] (5-66) 

Xl"'^ = 12y-^[{C + Dy) cosy + {Cy ~ D) siny - C] (5-67) 
A-f^ = j^y-'^[{2RC-Cy'^ + 2RDy) cosy + {2RCy-2RD + Dy'^) siny -2RC]i5-68) 

where y = x/y/S, and C and D are constants. There are also three solutions with vanishing 
Weyl tensor ($fc = 0), which we write as 

Zk = fA,y'''[2 + y^) (5-69) 
dfc = fA,y-' (5-70) 
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= -j^{AiCosy-A2smy + A^y-^) (5-71) 

= ^{-Aicosy + A,smy)-j=^Asy-' (5-72) 

X^'^^ = Aisiny + A2Cosy + A^y^'^ (5-73) 

= ^{A,smy + A2Cosy) + Asy-\ (5-74) 

where Ai, A2 and A3 are further constants. The solution with only ^3 non-zero describes 
a radiation dominated universe which is exactly FRW except that the CDM has a peculiar 
velocity (relative to the velocity of the fundamental observers) v^"^ = v^^ / S, where vj^^ is a 
first-order vector, orthogonal to the fundamental velocity, which is parallel transported along 
the fundamental flow lines: v^^ = 0(2). This can be seen most clearly by adopting the energy 
frame, defined by the condition = 0. This is arguably a better choice to make in the early 
universe, since u"" is then defined in terms of the dominant matter components, rather than a 
minority component, such as the CDM, which has little effect on the gravitational dynamics. 
Choosing the energy frame, and ignoring anisotropic stresses (which are frame-independent 
in linear theory), the CDM relative velocity evolves according to 



in the radiation dominated era. Since the CDM interacts with the other matter components 
through gravity alone, and since the gravitational influence of the CDM on the dominant 
matter components may be ignored during radiation domination, equation (|5-75| ) is the only 
equation governing the evolution of the perturbations which makes reference to the CDM. 
It follows that any solution of equation ( ^-75| ), defines a valid solution to the linearised 
perturbation equations. The solution, which has v^"^ = vj^^ /S" in an otherwise FRW universe 
corresponds to the solution labelled by A^, in equations ( |5-69| )-( [5-74| ). Note that this solution 



decays in an expanding universe. Following standard practice, we assume that this mode 
may be ignored (^3 = 0) since it would require highly asymmetric initial conditions at the 
end of the inflationary epoch if the decaying mode was significant during the epoch of interest 
here. Similar comments apply to the mode labelled by D in equations ( ^-62|) - (|5-6^ ). 

An important subclass of these solutions describe adiabatic modes. We assume that the 
appropriate covariant and gauge-invariant definition of adiabaticity is that 

for all species i and j (Bruni et al. 1992). This condition, which is frame-independent in 
linear theory, is the natural covariant generalisation of the (gauge- invariant) condition 

5p(*) 5p(^) 



p(i) + p{i) pU) + pU) ' 



(5-77) 
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where 5p 



P 



(i) _ nii) 



p^^' is the usual gauge-dependent density perturbation, and overbars 



denote the background quantity. Demanding adiabaticity between the photons and the 
neutrinos leaves only one free constant of integration, which we take to be C. The remaining 
constants are Ai = = D = and A2 = —6C. Note that the constants of integration will 
depend on the mode label k, in general, so we have C = Ck- 

This adiabatic solution may be developed further by including higher-moments of the 
neutrino distribution function, and finding a series expansion of the (enlarged) system in 



X. 



To obtain solutions correct to 0{x^), it is necessary to retain ir'j^'' and G)^' . The series 
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solution that results is 
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Note that on large scales (x ^ 1), the harmonic coefficient of Sab is constant along the 



flow lines. It follows that on large scales, we may write Sab 



SS!/S\ where S 



(0) 
ab 



0(2). 

The series solution given in equation ( p-7^ ) is adiabatic between the photons and neutrinos 
to O(x^), but the adiabaticity is broken by the higher-order terms. This difference in the 
dynamic behaviour of radiation and neutrinos is due to their different kinetic equations; the 
neutrinos are collisionless which allows higher-order angular moments in the distribution 
function to grow, but the radiation is tightly-coupled to the baryon fluid which prevents 
the growth of higher-order moments. The baryon relative velocity v^^ is determined by the 
condition that the radiation be nearly isotropic in the rest frame of the baryons: 



3 (7) 



(5-79) 



and the spatial gradients of the baryon and CDM follow from the adiabaticity condition: 



(7) 



(^) 



(5-80) 



where we have neglected the small effect of baryon pressure. The series solution given as 
equation ( [5-7^ ) was used to provide adiabatic initial conditions for the numerical solution of 
the perturbation equations, discussed in the next section. 



6. Adiabatic Scalar Perturbations in a = Universe 

In this section, we discuss the calculation of the CMB power spectrum from initially 
adiabatic scalar perturbations in an almost-FRW universe with negligible spatial curvature. 
The evolution of anisotropy in the CMB and inhomogeneities in the density fields, resulting 
from scalar perturbations, may be found by solving numerically the equations presented 
in Section |^, with initial conditions determined from the analytic solutions of the previous 
section. For adiabatic perturbations, the specification of initial conditions is particularly 
simple; there is a single function of the mode label k to set. This function gives the 
(constant) amplitude of the harmonic component of the electric part of the Weyl tensor on 
super-horizon scales. 



6.1. The CMB Power Spectrum 

The gauge-invariant temperature perturbation from the mean, denoted by StIg), is given 
by equation (|4-3|) . Substituting for the harmonic expansion of the angular moments Jal^^ai^ 
we find 



1 - (2/ + l)(20! 
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where we have introduced the radiation transfer function T^''\k), which is a function of the 
eigenvalue k only. The transfer function is defined to be the value of for the initial 
condition Ck = 1. Since the dynamics of a scalar mode labelled by the index k depends 
only on the eigenvalue of the eigenfunction Q^''\ the transfer function is a function of the 
eigenvalue k only. For the linearised theory considered here, we have = CkT'^^\k). 

We have come as far as we can without making a specific choice for the scalar harmonic 
functions Q*^*^-*. To proceed, we introduce an almost-FRW coordinate system [(Ellis 1996)| as 



follows. If the perturbations in the universe are only of scalar type, then the velocity 
is hypersurface orthogonal, so that we may label the orthogonal hypersurfaces with a time 
label t. Furthermore, since we have chosen m° to be the CDM velocity, which is geodesic, 
the flow orthogonal hypersurfaces may be labelled unambiguously with proper time along 
the flow lines, so that u"" = V°t. The orthogonal hypersurfaces depart from being spaces 
of constant curvature only at first-order, so we can introduce comoving spatial coordinates 
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x^, in such a way that our (synchronous) coordinate system is almost-FRW in form. (Latin 
indices, such as i, run from 1 to 3.) It is then straightforward to show that the functions 
^ik-x ^ where k x = and /c* are constants, satisfy the defining equations for the scalar 
harmonic functions, equations (5^) and ( p-7] ), with k"^ = k^k\ in an almost FRW universe 



with negligible spatial curvature. It follows that we may take 

Q(*=) = e**^-^. (6-2) 



For the open and fiat cases, the appropriate generalisations of the e*^'^ |^(Harrison 1967) 



should be used for the Q^''\ Note that the expansion coefficients, such as , depend on the 
detailed choice of the scalar harmonics Q^''\ but that covariant tensors, such as J2k Jk^Q^ai...ap 
are independent of this choice. If vector perturbations are also significant, we cannot use 
the velocity to define a time coordinate in the manner described above. Instead, an 
almost-FRW coordinate system should be constructed using an irrotational and geodesic 
velocity field m", which is close to our chosen fundamental velocity m". (One possibility is 
to take oc V"p.) Using this velocity field, almost-FRW coordinates can be constructed 
by the above procedure (Ellis 1996] . The resulting Q^^^ will satisfy the defining (zero-order) 



properties of the scalar harmonics in the frame, since the relative velocity of is first- 
order. 

Since u°'ea = 0, and = VH, we can always choose the so that at our observation 
point e'^ = (0, S'~^e*), with e*e* = 1 (for example, one can choose the so that Sx^ are 
locally Cartesian coordinates in the constant time hypersurface) . Then it follows that, to 
zero-order, 

Qil.,e-...e- = ^^|||£pK/^)Q^'=\ (6-3) 

where n = k'^e^/k, and Pi{fi) are the Legendre polynomials. This result demonstrates that 
expanding the angular moments of the distribution function in the covariant tensors Q^^ is 
equivalent to the usual Legendre expansion of the Fourier modes of the distribution function 
(which are axisymmetric about the wave vector k), in an almost-FRW universe, where the 
spatial curvature may be neglected. 

Following standard practice, we make the assumption that we inhabit one realisation of a 
stochastic ensemble of universes, so that the Ck are random variables. (The physical basis on 
which this assumption rests is that initial fiuctuations were generated from causal quantum 
processes in the early universe, such as during a period of infiation; see for example, Kolb 
& Turner (1990).) Given our chosen form for the Q^''\ statistical isotropy of the ensemble 
demands that the covariance matrix for the Ck takes the following form: 



(6-4) 
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where C^ik) is the primordial power spectrum which is a function of the eigenvalue k. 
The 6kk' appearing in equation ( |6-4| ) is defined by J2k ^kk'^k = ^k', where Ak is an arbitrary 
function of the mode label k. The CMB power spectrum Ci is defined by Ci = {\aim\'^)i where 
aim are the coefficients in the spherical harmonic expansion of the temperature anisotropy 
(see Section ||). Substituting the harmonic expansion of the J'i} ai i^^o equation ([4-5|), and 
using the zero-order result 

^a,., _ _ _ /.«'^'Q(^) ^^Q«*„^ = ^^|)||^, (6-5) 

which follows from Q^^^ = e"^^'^, we find the familiar expression for the CMB power spectrum 
in terms of the transfer functions and the primordial power: 

Ci = 7r^ Jd Ink C\k) [t(') {k)\\ (6-6) 

We make the standard assumption that on large scales, the primordial power spectrum 
may be approximated by a power-law of the form C^ik) oc k^"~^. Many inflationary models 
predict that the scalar index will be close to unity |(Kolb fc Turner 1990) . The case = 1 



describes the scale-invariant spectrum. This term arises from considering the logarithmic 
power spectrum in Fourier space of the (gauge-dependent) fractional density perturbation 
6p evaluated at horizon crossing. An analogous result can be found in the covariant and 
gauge-invariant approach: we evaluate the logarithmic power spectrum of the dimensionless 
scalar A = '^^)V^p/(pi^^) in the energy frame (g^ = 0). Making use of equation ( |2-15| ), with 



the contribution from anisotropic stress neglected, and the frame-invariance of ^^^'V"'Sab in 
linear theory, we flnd that 

d,^k{A') = fTrHj^'C'ik). (6-7) 

Note that A only receives a contribution at linear-order from scalar modes (Ellis et al. 1990). 
In deriving equation ( |6-71 ), which is valid before the modes labelled by k reenter the Hubble 
radius, we have assumed that only the fastest growing scalar mode is signiflcant so that $fc is 
constant before horizon crossing. For given k, the logarithmic power in A evolves in time due 
to the presence of Tik on the right-hand side of equation ( |6-7|) . However, at horizon crossing 
Tik falls below some critical value of order unity which is independent of k. It follows that 
for the scalar index = 1, the logarithmic power in A at horizon crossing is independent 
of scale. 



6.2. The Tight-Coupling Approximation 



At early times, when the baryons and photons are tightly-coupled, the radiation is 
nearly isotropic in the frame of the baryons. In this limit, it is convenient to replace the 
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propagation equations for the J^'-* with / > 1, and for the baryon relative velocity v^Ji'\ with 
approximate equations which may be developed by an expansion in the photon mean free 
path l/rieCTT- The approximate equations are simpler to solve numerically than the exact 
equations since the former do not include the large Thomson scattering terms present in the 
latter. 

For scalar perturbations, it is simplest to work directly with the harmonic expansion 
coefficients, and vl^\ The relevant time-scales in the problem are the photon mean free 
time tc = {ne(JT)~^i the expansion time-scale tn = H~^, and the light travel time across the 
wavelength of the mode under consideration, tk = S/k. In the tight-coupling approximation, 
we expand in the small dimensionless numbers tc/tn and tc/tk, so that the procedure is valid 
for tc <^ mm(tH,tk)- While a mode is outside the horizon, mm{tH,tk) = tn, whereas 
mm{tH,tk) = tk during the acoustic oscillations. In the CDM frame {wa = 0), the procedure 
is similar to that usually employed (Peebles & Yu 1970; Ma & Bertschinger 1995) in the 
synchronous- gauge. We combine the propagation equations given in Section ^]l] for the 
photon moments J^'"* {I > 1) and the baryon relative velocity vf'\ to get the exact (in linear 
theory) equations: 
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where = ql^^ — 4t>[^^/3, and for this section, R = Ap^'^y3p^^\ Iterating these equations 
gives the tight-coupling expansions 
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where the subscript on the variables on the right-hand side denotes the order in the expansion 
parameter e = ma.x{tc/tH,tc/tk). To avoid cluttering of indices, we leave the mode label k 
implicit. We shall only require the results to first-order in e: 
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The propagation equation for g^^^ in the tight-couphng approximation may be obtained from 



the exact equation (|6-8|) , with Ak replaced by ql 
their tight-couphng expansions. 



,(7) 



Av^^'' /3, and v^'' and n)^' replaced by 



.(7) 



6.3. Numerical Results 

We are now in a position to evolve an initial set of perturbations from early times to the 
present in an almost-FRW universe with negligible spatial curvature, and to calculate the 
power spectrum of the CMB anisotropies that result. In this section we present the results of 
a numerical simulation in the standard CDM model, with parameters Hq = 50kms~^Mpc~^, 
baryon fraction Qh = 0.05, CDM fraction = 0.95, Helium fraction 0.24, zero cosmological 
constant, and a scale-invariant spectrum of initially adiabatic conditions (n^ = 1). 

Our code to solve the covariant and gauge-invariant perturbation equations in the CDM 
frame, including the Boltzmann hierarchies for the photons and neutrinos, was based on 
the serial COSMICS code developed by Bertschinger and Bode, and described in Ma & 
Bertschinger (1995). The COSMICS package, including full documentation, is available at 
http://arcturus.mit.edu/cosmics. We modified the COSMICS code to solve the covariant 
equations given in this paper, for the matter variables and for the spatial gradient of the ex- 
pansion, Za- The shear, which is required to solve the Boltzmann hierarchies for the photons 
and neutrinos, was determined from the equation ( |2-lt: ). The electric part of the Weyl ten- 



sor. Sab, could then be determined from equation (|2-15|) . Our calculations of the zero-order 



ionisation history of the universe, which fully include the effects of Helium and Hydrogen 
recombination, followed Ma & Bertschinger (1995), as did our truncation schemes for the 
photon and neutrino Boltzmann hierarchies. The first-order tight-coupling approximation 
was used at sufficiently early times that max(tc/tH,tc/tk) <^ 1. 

In Figure ^ we show the variation of the harmonic coefficients X^^^ of the comoving 
fractional spatial gradients in the CDM frame, against redshift in the standard CDM model. 
Similar plots were given by Ma & Bertschinger (1995) for the Fourier components of the 
(gauge-dependent) density perturbations 5^*^ = (p(*) — p(*))/p(*), where p*^*^ is the density of 
the species i in the background model. Our results, given in Figure |l], agree well with the 
synchronous-gauge results of Ma & Bertschinger (1995). This is because the constant time 
surfaces in this gauge are orthogonal to the CDM velocity, so that A"^^*^ is a covariant measure 
of the density inhomogeneity in these surfaces. Although 5^*^ jg gauge-dependent in the 
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Fig. 1. — The variation of the harmonic coefficients of the fractional comoving density 
gradients in the CDM frame with redshift in the standard CDM model: — 0-95, = 0.05, 
Qa = 0, //o = 50kms-^Mpc-\ with Helium fraction 0.24, for k = 0.01, 0.1 and 1.0 Mpc"^ 
The normahsation is chosen so that = 1 for all k at early times. 
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synchronous gauge, the gauge-conditions restrict this gauge-dependence to transformations 
of the form 5^ _ ap/p, where a is a first-order constant. It follows that the Fourier 

coefficients of 5^ are gauge-invariant away from = in Fourier space. 

The qualitative behaviour of the comoving density gradients can be seen directly from 
their propagation equations. For scalar perturbations, it is simplest to work directly with 
the equations of motion ( |5-27| )-( ^-30| ) for the harmonic coefficients X^^^ in the CDM frame. 
Eliminating the spatial gradients of the expansion, Za, we find the following second-order 
equations: 

i 

+ 11^4") - nea^ (14^) - gi-)) (6-20) 

+ lext^ + l^xt^ = |. E {p^'^ + 3p^^0 - il^^i^^ + iS-i^^ (6-21) 

i 
i 

P^t^ + le^i:^ = IkJ:{p^'^+3p^^)4\ (6-23) 

i 

where we have ignored baryon pressure except in the acoustic term c1{k'^ / S'^)xl:^\ which can 
be significant on small scales. In the limiting case that the mode is well outside the Hubble 
radius (Tik ^ 1), the equations of motion for each component reduce to the common form 

A>« + |^#) = 4 (l + ?S) E {P^'^ + 4'^- (6-24) 

j 

For adiabatic initial conditions, it is clear that the adiabatic condition, given as equa- 
tion ( p-76| ), is maintained while the mode is outside the Hubble radius. Solving equa- 
tion ( 5-24| ) for adiabatic perturbations gives growing modes proportional to t and t"^^^ during 
radiation and matter domination respectively. 

If a mode enters the Hubble radius prior to last scattering, the photon/baryon fluid, 
which is still tightly-coupled, undergoes acoustic oscillations. To lowest-order in the tight- 
coupling parameter, ma.x{tc/tH,tc/tk), the photon and baryon perturbations remain adia- 
batic, evolving according to 

+3^^k +3{TTH}S^'^k -s^^l^yP J^fc - 9{T+:Rys^^fc ' 16-25) 

i 

where R = ip^^'^ /3p^''\ The solut ion of the homogeneous equation describes acoustic os- 
cillations in a fluid with sound speed squared {R + 3c^)/3(l -|- R), which are damped by 
the expansion of the universe. However, the oscillations are driven gravitationally by the 
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gradient *-^''Va(p + 3p), which gives an almost constant amphtude oscillation in the radia- 
tion dominated era. The Silk damping which is visible in Figure |I| for k = l.OMpc"^ at 
z ~ 10~^'^ arises from photon diffusion (which is not described by the lowest-order tight- 
coupling approximation) so is not described by equation ( 6-25 ). The neutrino perturbation 
also oscillates once inside the Hubble radius in the radiation dominated region, while the 
power-law growth of the CDM is impeded by the gravitational attraction of the oscillating 
dominant component (the inhomogeneous term in eq. [|6-23|| ). In the matter dominated era, 
the CDM becomes the dominant component, so we again see power-law growth of the CDM 
perturbation on all scales. Before last scattering, the photons and baryons remain tightly- 
coupled, but the character of the ^^^Va{p + 3]?) driving term in equation ( |6-25|) changes 
from an oscillation to a power-law as the CDM becomes dominant. At last scattering, the 
photons and baryons decouple. The baryons no longer feel the pressure support provided 
by the photons; the Jeans' length of the baryons is very small and the acoustic term in 
equation ( |6-22| ) is negligible. The ^^■'Va(p + 3p) driving term attracts the baryons into the 
potential wells caused principally by inhomogeneity of the CDM, so that Xj;^^ relaxes to Xj^'^^ 
as a power-law. After last scattering, the photons and neutrinos continue to undergo driven 
oscillations, which decay towards the particular integral X^'^^ = X^'^^ = QTi^xlf' . 

In Figure ^ we show the CMB power spectrum calculated from a simulation in the stan- 
dard CDM model. On large scales, the plateau arises from the usual potential fluctuations 
Y.k ^kQ^'^y^i on the last scattering surface [Sachs fc Wolfe 1967)| . The oscillations in the CMB 
power spectrum on smaller scales (the Doppler peaks) arise from the acoustic oscillations in 
the baryon/photon fluid. These oscillations give rise to strongly scale-dependent gradients 
of the photon energy density in the energy-frame, which in the approximation of instanta- 
neous recombination can be interpreted as temperature variations across the last scattering 
surface, and a local scale-dependent distortion of the last scattering surface relative to the 
energy frame. Since the last scattering surface is well approximated by a hypersurface of 
constant radiation temperature (so that recombination does occur there), it is more correct 
to interpret the Doppler peaks in terms of the local variations in redshift along null geodesies 
back to the last scattering surface, than in terms of temperature variations on the last scat- 
tering surface. (There is another significant contribution to the Doppler peaks, which is of 
dipole nature on the last scattering surface, and tends to fill in the power spectrum near 
the first Doppler peak; see Hu & Sugiyama (1995) and Challinor & Lasenby (1998) for more 
details.) On the smallest scales, the power spectrum is damped due to photon diffusion in 
the photon/baryon plasma prior to recombination. 
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Fig. 2. — The power spectrum of scalar CMB anisotropies in the standard CDM modeh 
fic = 0.95, fife = 0.05, = 0,Ho = 50kms~^Mpc-^ n, = 1, with Hehum fraction 0.24. The 
normahsation of the vertical scale is arbitrary. 

7. Tensor Perturbations 

The covariant equations of Section ^ are independent of any non-local splitting of the 
perturbations into scalar, vector and tensor modes. In the linear approximation, we have 
seen that the equations describing scalar perturbations can be obtained from the full set 
of equations in a straightforward manner. The same is also true for vector and tensor 
perturbations. Vector perturbations are not expected to be important today in non-seeded 
models since the vorticity decays in an expanding universe [(Hawking 1966) . However, most 



inflationary models do predict the generation of a primordial spectrum of gravitational waves 
(tensor modes) during an epoch of inflation (see Lidsey et al. (1997) for a comprehensive 
review). The existence of such a background of relic gravitational waves would have a 
significant effect on the CMB anisotropy power spectrum on large scales [(Crittenden et aL 



1993 )[ . For completeness, we give the tensor multipole equations in this section. We defer 
a detailed derivation of the equations, and a discussion of their relation to the multipole 
equations usually employed in calculations of the effects of tensor modes [(Crittenden et al.] 



1993 )[ to a future paper [(Challinor 1998^ 
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In the covariant approach to cosmology, gravitational waves are characterised by re- 
quiring that the vorticity and all gauge-invariant vectors and scalars vanish at first-order 
(Dunsby et al. 1997), so that the spatial gradients of the density and expansion vanish, as 
well as the acceleration and the heat fluxes. The electric and magnetic parts of the Weyl 
tensor, the shear tensor and the anisotropic stress are constrained to be transverse: 

(3)V"£,, = 0, (^V'^S,, = 0, VV,, = 0, (=^VV,, = 0. (7-1) 

It is straightforward to show that these conditions are consistent with the linearised propa- 
gation equations given in Section ^ (see also Maartens (1997) for a discussion of consistency 
of the exact equations for irrotational dust spacetimes). As with scalar perturbations, it 
is convenient to expand harmonically the first-order covariant variables in tensors derived 
from solutions of a generalised Helmholtz equation. For tensor perturbations, we employ the 



tensor- valued solutions [(Hawking 1966)| Qit^ = Q/qM (note that we use the same symbol for 



the tensor harmonics and the second-rank tensors derived from the scalar harmonic functions 
to avoid cluttering formulae with additional labels), which satisfy the zero-order relations 

gltV = 0, (3)vag(J)^Q^ Q1?=0. (7-2) 

Expanding Sab, ^ab, o'ab and Tiab in these tensors (as in Section | but with Q^^'' replaced by 
the tensor harmonics, and the scalar- valued variables $fc replaced by Sk), we obtain simple 
propagation equations for the electric part of the Weyl tensor and the shear: 

+ + f -i7/^p)cTfc = -l(37-l)/tp^7r, + i/€p7T, (7-3) 

|((Tfe + i^(Tfe) + (|)'^fc = -Inpn,. (7-4) 



Note that we have used the constraint equation ( p-13| ) to eliminate the magnetic part of the 



Weyl tensor from the propagation equation for the electric part. Equations ( |7-3[ ) and ( [7-4| ) 
close up, with the anisotropic stress treated as a known field, to give a second-order equation 
for the shear: 

ak + 0&k+[f + ^- 1(37 - 2)kp\ ak = np^ [|(37 - 2)07r,. - tt^] , (7-5) 

which generalises the homogeneous equation derived in Dunsby et al. (1997) to include 
anisotropic stress. 

For the photon and neutrino angular moments, J'^^ and G^'^-* , we expand in tensors 
Q^ai ai derived from the tensor harmonics using the same recursion relation as for scalar 
perturbations, 

r){k) _ 5 A3)y _ j^(3)yf>^(fc) , \ 

Vai...a, - fc 1^ ^(aiVa2...a,) 2l~l ^ '^b(ai...ai_2'^ai_iai) ) , I'-Oj 
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for / > 2. This procedure gives the following covariant Boltzmann multipole equations for 
tensor perturbations: for the anisotropic stress (/ = 2), 

+ 'si - W) 4'' - M-^ = ->ea.vrr (7-7) 
+ -^l-^' = 0, (7-8) 



and for I > 3, 

4'' + l{|Si5%[i-((' + i)'-3)fl4'""-^4'-"} = -".-t4" (7-9) 

'^k + 5 I (Z+l)(2«+l) ~ + -'^jFj'-^fc ~ 2m'-^fc i ~ 

As we stressed earlier, by performing the covariant angular expansion before harmonically 
expanding the moments, the necessary angular dependence of the moments appears auto- 
matically, whereas in the metric-based approach, this (azimuthal) dependence must be put 
in by hand, and is different for the two polarisations of gravitational waves. It should be 
noted that the moment equations given here are not the same as those satisfied by the A}; 
variables, where i = +, x labels the polarisation of the gravitational wave and / denotes the 
intensity Stokes parameter (we follow the notation employed in Kosowsky (1996)), that are 
usually employed in metric-based calculations of the effects of tensor modes. This is because 
the covariant angular expansion gives rise to a more natural set of variables, the , which 
are related to the temperature anisotropy in a simpler manner than the A};. In particular, 
the l-th multipole Ci of the anisotropy power spectrum depends only on whereas Ci 
depends on the (Z — 2), / and (/ + 2)-th moments. A};, thus obscuring the physical interpre- 
tation of these variables. The relation between the two sets of variables will be discussed 
further in Challinor (1998). 

The first-order propagation equations for the shear and electric part of the Weyl tensor, 
along with the Boltzmann multipole equations for the photons and neutrinos give a closed 
set of equations that can be solved to calculate the temperature anisotropy for given initial 
conditions. The numerical solution of these equations will be considered elsewhere |( Challinor 
119981. 



8. Conclusion 

We have shown how the full kinetic-theory calculation of the evolution of CMB anisotropics 
and density inhomogeneities can be performed in the covariant approach to cosmology 
(Ehlers 1993; Ellis 1971), using the gauge-invariant perturbation theory of Ellis & Bruni 
(1989). Adopting covariantly-defined, gauge-invariant variables throughout ensured that 
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our discussion avoided the gauge ambiguities that appear in certain gauges, and that all 
variables had a clear, physical interpretation. We presented a unified set of equations de- 
scribing the evolution of photon and neutrino anisotropics and cosmological perturbations 
in the CDM model, which were independent of a decomposition into scalar, vector or tensor 
modes and the associated harmonic analysis. Although we only considered the case of linear 
perturbations around an FRW universe here, it is straightforward to extend the approach to 
include non-hnear effects (Maartens et al. 1998), which should allow a physically transparent 
discussion of second-order effects on the CMB. Indeed, the ease with which one can write 
down the exact equations for the physically relevant variables is one of the major strengths 
of the covariant approach. 

The linear equations describing scalar modes and tensor modes were obtained from the 
full set of equations in a straightforward and unified manner, highlighting the advantage of 
having the full equations (independent of the decomposition into scalar, vector and tensor 
modes) available. For the scalar case, the Boltzmann multipole equations for the moments 
of the distribution functions obtained here were equivalent to those usually seen in the liter- 
ature. However, for tensor modes, the covariant approach led naturally to a set of moment 
variables that more conveniently describe the temperature anisotropy than those usually 
employed. For scalar modes, we discussed the solution of the perturbation equations in 
detail, including the integral solution of the Boltzmann multipole equations and the rela- 
tion between the timelike integrations performed in the multipole approach to calculating 
CMB anisotropies, and the lightlike integrations of the fine of sight approach. The numer- 
ical solution of the scalar equations in a X = 0, almost-FRW, CDM universe were also 
discussed. Our numerical results provide independent confirmation of those of other groups, 
(see, for example. Ma & Bertschinger (1995) and Seljak & Zaldarriaga (1996)), who have 
obtained their results by employing non-covariant methods in specific gauges. Typically, 
these methods require one to keep careful track of all residual gauge-freedom, both to enable 
identification of any gauge-mode solutions, and to ensure that the final results quoted are 
gauge-invariant (and hence observable). Fortunately, the isotropy of the photon distribution 
function in an exact FRW universe ensures that the CMB power spectrum, as calculated 
from the gauge-dependent perturbation to the distribution function, is gauge-invariant for 
I > 1. 

We hope to have shown the ease with which the covariant approach to cosmology can be 
applied to the problem of calculating CMB anisotropies. The covariant and gauge-invariant 
method discussed here frees one from the gauge problems that have caused confusion in 
the past, and focuses attention on the physically relevant variables in the problem and the 
underlying physics. Future work in this area will include the discussion of non-linear effects 
(Maartens at al. 1998), the inclusion of polarisation, and the effects of hot dark matter, all 
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of which can be expected to bring the same advantages of physical clarity and transparency 
that we hope to have demonstrated here. 

The development of the COSMICS package was supported by the NSF under grant 
AST-9318185. The authors wish to thank Roy Maartens for useful comments on an earher 
version of this paper. 
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